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Abstract 

Aiming at the study of critical phenomena in the presence of boundaries with a non-trivial shape 
we discuss how lattices with an adaptive lattice spacing can be implemented. Since the parameters 
of the Hamiltonian transform non-trivially under changes of the length-scale, adapting the lattice 
spacing is much more difficult than in the case of the numerical solution of partial differential 
equations, where this method is common practice. 

Here we shall focus on the universality class of the three-dimensional Ising model. Our starting 
point is the improved Blume-Capel model on the simple cubic lattice. In our approach, the system 
is composed of sectors with lattice spacing a, 2a, 4a, ... . We work out how parts of the lattice with 
lattice spacing a and 2a, respectively, can be coupled in a consistent way. Here, we restrict ourself 
to the case, where the boundary between the sectors is perpendicular to one of the lattice-axis. 
Based on the theory of defect planes one expects that it is sufficient to tune the coupling between 
these two regions. To this end we perform a finite size scaling study. However first numerical 
results show that slowly decaying corrections remain. It turns out that these corrections can be 
removed by adjusting the strength of the couplings within the boundary layers. 

As benchmark, we simulate films with strongly symmetry breaking boundary conditions. We de¬ 
termine the magnetization profile and the thermodynamic Casimir force. For our largest thickness 
Lq = 64.5, we find that results obtained for the homogeneous system are nicely reproduced. 
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I. INTRODUCTION 


The present work is motivated by the numerical study of effects related to critical phenom¬ 
ena in the presence of boundaries with non-trivial shapes. Our own interest in this problem 
stems from an ongoing project on the thermodynamic Casimir effect involving hnite 
objects. In particular, in ref. we studied the thermodynamic Casimir force between a 
sphere and a plane by using Monte Carlo simulations of the improved Blume-Capel model 
on a simple cubic lattice. 

In the numerical study of partial differential equations often adaptive grids are used. In 
regions, where the held varies strongly, a hner resolution of the space is used than in regions, 
where the held is smooth. In fact, solving problems related to the thermodynamic Casimir 
ehect in the mean-held approximation, a hnite element method with a mesh adapted to the 
problem is used. For recent work on three-body interactions see for example refs. @,0- 

In the case of critical phenomena, beyond the mean-held approximation, this approach 
is complicated by the fact that helds and the strength of the interaction transform in a 
non-trivial way under a rescaling of the length. Therefore it is hard to solve the problem in 
general. Here we consider the Ising universality class in three dimensions and in particular 
study the improved Blume-Capel model. We have chosen a rather simple setup. We start 
from a simple cubic lattice with the lattice spacing a. Certain regions of the lattice are coarse 
grained and the lattice spacing is increased to 2a, 4a, ... . These regions are composed of 
cubes with faces parallel to the (100), (010), and (001) planes. Here we perform a hrst step 
of this program. We study how to couple two half spaces with lattice spacing a and 2a, 
respectively. To this end, we employ finite size scaling to determine the proper values of the 
couplings at the boundary between the two different lattice spacings. 

In principle the boundary between lattice spacing a and 2a can be seen as a defect plane. 
As argued in refs. 0 0 such a defect is a relevant perturbation of the non-trivial hxed point 
that characterizes the continuous phase transition of the bulk systems. Its RG-exponent is 
Ud = Ut ~ ^1 where yt = 'Vjv and v is the critical exponent of the correlation length of the 
bulk system. This suggests that it is sufficient to tune the coupling between the regions 
with lattice spacing a and 2a. In our hnite-size scaling study we consider a lattice of the 
size Lq X Li X L 2 , where Lq = 2L and Li = L 2 = L. We split the lattice in two halfs in 
0-direction. In one half we use lattice spacing a, while in the other one the lattice spacing is 
2a. We define phenomenological couplings that are suitable to detect correlations between 
the two parts of the lattice with different lattice spacing. First numerical results show that 
there are slowly decaying corrections. It turns out that these can be eliminated by tuning 
the couplings within the lattice planes at the boundary. 

Finally we benchmark our result by simulating films with strongly symmetry breaking 
boundary conditions. We replace half of the sites, located in the middle of the film, by 
coarse grained ones. We compute the thermodynamic Casimir force and the magnetization 
profile. These are compared with estimates that we obtained in ref. for the homogeneous 
system. 

The paper is organized as follows. First we define the Blume-Capel model and summarize 
relevant previous work. Next we introduce the finite size scaling method used to determine 
the couplings at the boundary between different lattice spacings. Then we discuss our 
simulations and the numerical results. We summarize our estimates for the couplings at the 
boundary between different lattice spacings. Then we present our simulations of dims with 
strongly symmetry breaking boundary conditions. Finally we conclude and sketch how we 
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shall use our outcome in upcoming studies of the thermodynamic Casimir force. 


II. THE MODEL 

We consider the Blume-Capel model on a simple cubic lattice. The reduced Hamiltonian 
is given by 

H=-/3J2 + D^sl-h^s^ , ( 1 ) 

<xy> X X 

where Sx G {—1, 0,1} and x = (xq, xi,X 2 ) denotes a site on the simple cubic lattice, where Xi 
takes integer values and < xy > is a pair of nearest neighbors on the lattice. The partition 
function is given by Z = where the sum runs over all spin conhgurations. 

The parameter D controls the density of vacancies s^: = 0. In the limit D ^ —oo vacancies 
are completely suppressed and hence the spin-1/2 Ising model is recovered. Here we consider 
a vanishing external held h = 0 throughout. In d > 2 dimensions the model undergoes a 
continuous phase transition for —oo < D < Dtri at a f3c that depends on D. For D > Dt^i 
the model under goes a hrst order phase transition, where for the three-dimensional lattice 
Dtr^ = 2.0313(4) ]n|. 

Numerically, using Monte Carlo simulations it has been shown that there is a point 
{D*, PdD*)) on the line of second order phase transitions, where the amplitude of leading 
corrections to scaling vanishes. Our recent estimate is D* = 0.656(20) [l^. In [l^ we 
simulated the model at D = 0.655 close to Pc on lattices of a linear size up to L = 360. From 
a standard hnite size scaling analysis of phenomenological couplings like the Binder cumulant 
we hnd /9c(0.655) = 0.387721735(25). Furthermore the amplitude of leading corrections to 
scaling is at least by a factor of 30 smaller than for the spin-1/2 Ising model. In the following 
study we shall simulate at D = 0.655 throughout. 


III. THE FINITE SIZE SCALING METHOD 

We study lattices of the size Lq = 2L and Li = L 2 = L with periodic boundary conditions 
in all three directions. In order to get numbers to compare with, we simulated a lattice with 
a unique lattice spacing a = 1 hrst. 

Next we studied a lattice with the lattice spacing a = 1 for 0 < Xq < T and lattice 
spacing 2a = 2 for L < xq < 2L. In the following we shall refer to the hne and the coarse 
part of the lattice simply by hne or coarse lattice, respectively. The sites of the hne lattice 
are located at xq = 0.5,1.5,..., L — 0.5, while those of the coarse lattice assume the positions 
Xq = L + 1,L + 3,...,2L — 1. Between the two parts of the lattice there are two boundaries 
that are perpendicular to the 0-axis and are located at xq = 0 and xq = L. 

The coupling constant in the bulk of the hne lattice is called Pi, the coupling constant 
in the bulk of the coarse lattice is called P 2 - The coupling between the hne and the coarse 
lattice is denoted by Pi^ 2 - It turned out that also the couplings within the layers next to 
the boundaries have to be tuned. They are denoted by Pi^b and P 2 ,b, respectively. We 
parameterize these couplings by Pi^b = Ci^bPi- This means that the reduced Hamiltonian of 
our system is given by 

H = —Pi SxSy — P2 SxSy — Pl,b SxSy — P 2 ,b SxSy 

<xy>eli <XY>el2 <xy>£Bi <XY>&B2 
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FIG. 1. We give a two-dimensional sketch of the lattice with two lattice spacings. The couplings 
in the bulk of the fine and the coarse lattice are denoted by j3i and /32, respectively. The couplings 
within the layers at the boundary between the two parts of the lattice are denoted by j3i^b and /32,fe, 
respectively. The corresponding links are given by dash-dotted lines in the sketch. The coupling 
between the hne and the coarse lattice are denoted by /3i,2- The corresponding links are indicated 
by blue dashed lines. 


-/ 3 i .2 + £> 5 ^ + £> 5 ^ si ( 2 ) 

<xY> X X 

where Ji and I 2 denote the set of links in the interior of the hne and the coarse lattice, 
respectively. Bi and B 2 denote the set of links within the layers next to the bonndaries. 
The snbscripts 1 and 2 indicate the hne and the coarse lattice, respectively. Sites of the hne 
lattice are denoted by x and ?/, while those of the coarse one are denoted by X and Y. A 
site in the coarse lattice, adjacent to the bonndary, is conpled to four sites of the hne one. 
The spins on the coarse lattice are Sx G In Fig. [T]we give a two-dimensional 

sketch. 

Next we dehne the quantities that we computed. We consider the magnetization in the 
center of the hne lattice 

Ml = Sx (3) 

L/4<xo<3L/4,xi,a:2 

and the magnetization in the center of the coarse lattice 

M2 = Sx ■ (4) 

L+L/A<Xo<L+'iL/A,Xi,X2 


It is plausible that the normalized correlation function 


(M 1 M 2 ) 

YMml) 


(5) 


is well suited to tune the coupling /3i^2 between the two lattices. In addition we determined 
the fourth order cumulant 


4". 


{Mil 

{Mf)^ 


( 6 ) 


4 






































where i G {1, 2}. 

We study the expectation value of the square of the magnetization 





( 7 ) 


where the normalization is chosen such that q{xQ) varies slowly with L. At the critical point 
we get 

g(xo) oc L"’' . (8) 

Furthermore we consider the prohles of the energy density 


Eo{xq) 


1 

1 ? 


E' 


Xl,tC2 


^XQ,Xl,X2 


'^XO + 1 


X\,X2l 


(9) 


and 


^ 2 ( 2 ^ 0 ) JJ2, ^ ^ {Sxq,x\,X2-^xq,x^+1,X2 a ^xo,xi,X2^xo,xi,X2+l) 


Xi,X2 


( 10 ) 


Analogous dehnitions are used for the coarse lattice. In the case of a homogeneous system, 
g(xo), Eo{xo), and Eui^xo) do not depend on xq. In the system with lattice spacings a = 1 
and 2a = 2 for the two parts of the lattice, this behavior should be reproduced as good as 
possible. Deviations should rapidly decay with the distance from the boundaries between 
different lattice spacings. 


IV. NUMERICAL RESULTS 


A. Simulations at the critical point 


First we simulated at our best estimate of the critical coupling f3i = (32 = 0.387721735, 
ref. 1^. As a preliminary step we determined the quantities 6 * 1 , 2 , V, q, Eq, and Eu dehned 
above for lattices with a unique lattice spacing a = 1 for L = 8 , 12, 16, 24, ..., and 128. For 
equilibration we performed 10 ® times two sweeps with the local heat-bath algorithm followed 
by L single cluster updates [l3| . Performing a sweep means that we run through the lattice 
in type-writer fashion. Then we evaluated the quantities of interest after two sweeps with 
the local heat-bath algorithm followed by L single cluster updates [l^ . In the following we 
shall refer to the evaluation of the quantities as performing a measurement. The number 
of these measurements for each lattice size is given in table [H where we also summarize our 
estimates of the observables. In total these simulations took about 8 years of CPU-time on 
a single core of a Quad-Core AMD Opteron(tm) 2378 CPU. 

Then for systems with two different lattice spacings we tuned the couplings at the bound¬ 
ary between the two lattice spacings such that we reproduce the values of the quantities ( 61,2 
and Vi- The quantities q, Eq and E 12 are studied to further validate the results obtained. 
Also here we performed a measurement after two sweeps with the heatbath algorithm fol¬ 
lowed by L single cluster-updates. As seed of the cluster we take with probability 1/3 a site 
of the coarse lattice and with probability 2/3 a site of the hne lattice. This way the average 
volume of the clusters in the hne and the coarse lattice is approximately equal. In addition 
to the quantities dehned in eqs. fl5]l6ll7ll9][T0|) we have implemented the correlators needed 
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TABLE I. Results obtained for lattices with a unique lattice spacing at (5 = 0.387721735. The 


value of the energy density in the limit L —>■ oo is taken from ref. "l4^. In the first column we give 


the lattice size L. In the second column we give the number of measurements divided by 10®. In 
the following columns we give our numerical estimates of the quantities defined in eqs. (1516171 iillOl) , 
respectively 


L stat/10® Gi ^2 

12 1005 0.621630(27) 

16 1015 0.623497(27) 

24 1008 0.624844(28) 

32 921 0.625336(30) 

48 816 0.625699(34) 

64 380 0.625721(50) 

96 302 0.625975(59) 

128 270 0.626055(64) 


V 

1.794209(46) 

1.794036(47) 

1.793880(50) 

1.793806(54) 

1.793627(60) 

1.793701(90) 

1.79355(11) 

1.79347(12) 


Eq 

0.2126616(18) 

0.2087201(13) 

0.20524821(85) 

0.20373756(64) 

0.20241808(43) 

0.20184545(45) 

0.20134865(31) 

0.20113357(23) 


_ Eu _ 

0.2129622(18) 

0.2088462(13) 

0.20528524(85) 

0.20375324(64) 

0.20242275(43) 

0.20184739(45) 

0.20134924(31) 

0.20113379(23) 


Q 

0.933959(32) 

0.926344(33) 

0.914493(36) 

0.905711(39) 

0.893115(43) 

0.883943(65) 

0.871432(76) 

0.862552(83) 


to compute the first derivative of these quantities with respect to j3i^2- With hindsight we 
would have also implemented the derivatives with respect to the other couplings. For given 
and C 2 , 6 , we determined the proper value of (3i^2 by requiring 

G,a(L, A,2) = G(L) := („) 


where the subscript horn indicates the homogeneous system. With our dehnition of G{L) 
we aim at a cancelation of corrections to scaling to a large extend. In the following we shall 
indicate results taken from the homogeneous system by a bar. We solved eq. flTT]) iteratively. 
Starting from a hrst guess jdfl for the solution of eq. flTT]) we iterated 


GigLjfl) 


+ 


dGi,, 

d/3i,2 


h+i) 

1,2 


dS )= g(l) 


( 12 ) 


until small. The number of measurements at a given is increased as we 

approach the solution. For larger lattice sizes we took from the results that we already 

obtained for smaller ones. 

First we simulated at = 02 ^ = 1 systems with the lattice sizes L = 32, 64 and 128. At 
the resulting values of ^ 1^2 we computed also the other quantities of interest. In Fig. [2] we 
plot q{xo)/q, where q is the estimate obtained for the homogeneous system with the lattice 
size L and L/2 for the hne and the coarse lattice, respectively. 

We see that there are large corrections at the boundaries. For the hne lattice q{xQ) is too 
small and too large for the coarse lattice. This effect only slowly decreases with increasing 
lattice size. This rather strong deviation from the behavior of the homogeneous system can 
also be seen clearly in the estimates of Vi and V 2 reported in table HIl We also give the 
estimates of 

flv = F - 1 (13) 

V2 
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FIG. 2. Profiles, eq. d?]), at ci^j, = C 2 ^b = 1 and the lattice sizes L = 32, 64 and 128. We plot 
q = q{xo)/q as a function of x = xq/L and x = {xq — L)/L for the fine and the coarse lattice, 
respectively. The curves for the different lattice sizes can be identified by the number of data 
points. Furthermore ^ < 1 for the fine lattice and q > 1 for the coarse one. The dashed lines 
should only guide the eye. 


TABLE II. Results for I3i^2, Vi smd V 2 at ci^b = C 2 ^b = 1 and /3c- In the last column we give Ry as 
defined in eq. m- 

~ h,2 Vi-V V2-V Ry 

32 0.191000(23) 0.0213(3) -0.0265(3) 0.0262(2) 

64 0.190598(48) 0.0169(8) -0.0194(7) 0.0200(6) 

128 0.190505(25) 0.0121(7) -0.0147(6) 0.0147(5) 


which allows us to monitor the asymmetry between the fine and the coarse lattice. 

We see that \Vi — V\, \V 2 — V\, and \Rv\ are decreasing rather slowly with increasing 
lattice size L. Fitting with the Ansatz 


Vi-V = cL-‘^^ (14) 

we obtain uy ~ 0.4. For V 2 — V and Ry we get similar results. 

It turns out that this deviation from the desired behavior can be lifted to a large extend 
by tuning ci^b and C 2 fi. We employ the quantities Vi and V 2 to this end. By construction, 
they are not sensitive to corrections that decay rapidly with increasing distance from the 
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TABLE III. Estimates of the optimal {ci^, C 2 ,b) at j3i^2 s-nd /3i = /32 = Pc as a function of the lattice 
size L. We also give the sum and the difference of the optimal ci^h and C 2 fi- 


L Cl^b _ C2,b Cl^b + C2,b Cl^b — C2,b 

24 1.0685(5) 0.9237(5) 1.9922(8) 0.1448(3) 

32 1.0679(6) 0.9233(7) 1.9912(11) 0.1446(3) 

48 1.0681(9) 0.9235(12) 1.9916(20) 0.1446(4) 


boundary between the fine and the coarse lattice. 

First we searched for choices of ci^b and C 2 ,b, where Ry vanishes. Our simulations show that 
this is the case for a line in the {ci^b, C 2 , 6 )-plane. At the level of the precision of our preliminary 
study this is the case for e.g. {ci^b,C 2 ,b) = (1,0.843), (1.072,0.93), (1.128,1). However, it 
turns out that in general still Vi ~ V 2 deviate from the value obtained for the homogeneous 
system. For (ci^t, C 2 ,b) = (1,0.843) we get Vi — V = 0.0052(3), 0.0025(4), and 0.0021(6) for 
L = 32, 64, and 128, respectively. On the other hand, for (ci^fe, C 2 ,b) = (1.128,1) we get 
Vi-V = -0.0059(3), -0.0028(4), and -0.001_9(5) for L = 32,’ 64,’and 128, respectively. 
For (ci^b,C 2 ,b) = (1.072,0.93) we find that Vi — V approximately vanishes. 

Motivated by this preliminary study, we systematically searched for the (ci,b, C 2 ,b), where 
both Vi — V and V 2 — V vanish at Pi^ 2 - To this end we simulated systems with L = 24, 
32 and 48 for various choices of (ci^t, 02 ^) and interpolated Vi — V and V 2 — V linearly in 
(ci,6, C 2 ^b)- Our results are summarized in table IIIII It turns out that the difference Ci^b — C 2 ,b 
has a much smaller error than the sum Ci^b + C 2 ,b. This might correspond to the fact that 
for the line where Ry = 0 the difference ci^b — C 2 ,b is roughly constant. 

Since the estimates obtained for L = 24, 32 and 48 are consistent within errors, we 
abstained from studying larger lattice sizes at this point. As our final choice we take 

(ci,b,C2,fe) = (1.068,0.9234) . (15) 

We had simulated at L = 24, 32 and 48 at {cyb, C 2 ,b) = (1-068, 0.923) and (1.068, 0.924) 
with a high statistics of 10® measurements or more. In table HV] we interpolate the numerical 
results to (ci^b, C 2 ,fe) = (1.068,0.9234). Furthermore we report the results of simulations per¬ 
formed directly at (ci^b, C 2 ,b) = (1-068, 0-9234) for L = 64, 96 and 128- For these lattice sizes 
we performed 94 x 10®, 61 x 10® and 46 x 10® measurements, respectively- The simulations 
of the L = 128 lattice took about 5-4 years of CPU-time on a single core of a Quad-Core 
AMD Opteron(tm) 2378 CPU- 

We find that for L = 64, 96 and 128 the estimates of Vi — U and V 2 — V remain consistent 
with zero- We fitted our data for Pi ^2 with the Ansatz 

^i, 2(T) = ;di^2(oo)-|-cL (16) 

Taking e = 2 and the data for L = 48, 64, 96 and 128 we arrive at /3i^2(oo) = 0-1858085(39), 
c = 0-054(14) and x^/d- 0 -f-= 0-36- In order to check the dependance of the result on e 
we repeated the fit with e = 1-6- We arrive at /3i^2(oo) = 0-1858061(45), c = 0-013(3) and 
X^/d- 0 -f-= 0-38- Since the estimates of Pi ^2 obtained for L = 64, 96 and 128 are consistent 
within errors, it is not surprising that /3i^2(oo) does not depend strongly on e. As our final 
estimate we take 


/?i,2(oo) = 0-185807(6) - 


( 17 ) 





TABLE IV. Results for /3i^2, Vi — V and V 2 — ^ at the optimal (ci^b,C 2 ,b) = (1-068,0.9234) at 
/3i = /32 = /3c- Ror a discussion see the text. 


L 


A,: 


Vi-V 


24 0.1858875(47) 
32 0.1858599(43) 
48 0.1858328(39) 
64 0.1858187(46) 
96 0.1858184(45) 
128 0.1858099(45) 


0 . 00010 ( 6 ) 

-0.00003(6) 

0.00001(7) 

-0.00002(11) 

-0.00007(13) 

-0.00001(15) 


V 2 -V 

-0.00000(5) 

0 . 00000 ( 6 ) 

0 . 00001 ( 6 ) 

0 . 00002 ( 8 ) 

0.00000(9) 

0.00002(13) 


where the error that we quote also takes into account the systematic uncertainty of the 
extrapolation to L —>■ 00 . 

In Fig. [3] we plot q = q{xo)/q, eq. ([7]), for (ci^b, C 2 ,fe) = (1.068,0.9234), /3 = 0.387721735 
and /3 i 2 = 0.185807. The improvement compared with (ci^b, C 2 ,fe) = (1,1) potted in Fig. [2]is 
fairly obvious. Starting from a distance of roughly 10 lattice spacings from the boundary we 
hnd that q{xo) is consistent, at the level of our statistical accuracy, with the result obtained 
for the homogeneous system. This observation holds for both the hne as well as for the coarse 
lattice. We checked also the behavior of Eq{xo) and Ei 2 {xo), eqs. fl^lITU]) . as a function xq. 
Here we hnd that starting from a distance of roughly 13 lattice spacings from the boundary 
Eq{xo) and Ei 2 {xo) are consistent with the results obtained for the homogeneous system. 

To be on the safe side, studying properties of the homogeneous system, one should keep 
at least a distance of ~ 20 lattice spacings from the boundary between different lattice 
spacings. 


B. Going ofF-critical 


The coupling /32 in the bulk of the coarse lattice as a function of the coupling /3i of the 
hne lattice is chosen such that the bulk correlation length is the same for both parts of the 
lattice. Hence in terms of the correlation length ^ in units of the respective lattice spacing 
we get 

e(/3i) = 2 ^/ 32 ) . (18) 

In the neighborhood of the critical point, the correlation length behaves as 

m = {i+b\tf+ct+d\tf+e\tr+...), ( 19 ) 


where t = (3c — (3 is our dehnition of the reduced temperature. Here we discuss an improved 
model which is characterized by 5 ~ 0. Hence, in the following we shall neglect the terms 


b\t\^ and Note that 6' = z/o;', where uj' = 1.67(11), ref. [15|. Hence in the numerical 


analysis of data for the correlation length it is difficult to disentangle the corrections ct and 
d\t\^ . After a preliminary study, we realized that for our purpose corrections to scaling can 
be ignored entirely. Hence we get the linear relation 


(ft - ft) = (ft - /3i) a'''" , 


( 20 ) 
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FIG. 3. Profiles, eq. dZ]), at ci,b = 1.068, C 2 ,b = 0.9234, /3 = 0.387721735, A ,2 = 0.185807 and the 
lattice sizes L = 64 and 128. We plot q = q{xo)/q as a function of x = xq/L and x = {xq — L)/L 
for the fine and the coarse lattice, respectively. The curves for the different lattice sizes can be 
identified by the number of data points. Furthermore q < 1 for the fine lattice and q > 1 for the 
coarse one. The dashed lines should only guide the eye. 


where we take u = 0.63002, ref. [12|. In order to keep things simple, ci^;, and C 2 ,b keep their 


valnes as determined above at the critical point. 

It remains to determine /3i,2 as a fnnction of (3i. Analogous to eq. flTT]) we impose 


Gi,2{L, f3i, f3i^2) — G{L,(3 i) 


Gl^2,hom{L/2, (32) + Gi^2,hom{.L^ (3i) 
2 


( 21 ) 


We have simulated L = 24 and 48 at one value of /5i in the high temperature phase and one 
value of /9i in the low temperature phase. In order to check the dependance of the result 
on the difference of these two values of /9i, we simulated L = 32 at four values of /9i. The 
corresponding results for the slopes A/3i^2/A/3i are given in table IVl 

Within the numerical accuracy we see no non-linearity nor a dependence on L. As our 
hnal result we take 

= 0.974(10) (22) 

dpi 


at the critical point. 


10 















TABLE V. Results for the slopes A/3i^2 /at (ci^;,,C 2 ,b) = (1.068,0.9234) 


L /3i,high /3i,low A/3 i,2/A/3i 
24 0.386 0.39 0.9738(38) 

32 0.386 0.3895 0.9706(51) 

32 0.387 0.3885 0.9774(68) 

48 0.3871 0.3883 0.9663(99) 


C. Summary of our results 


The relation of the bulk coupling on the hue f3i and the coarse lattice (32 is given by 


eq. 


(ft- A) = (A-ft) a*-'- , 


(23) 


where we take /3c = 0.387721735 and u = 0.63002, ref. [12|]. The amplification of the 
coupling in the layers at the boundary of the coarse and the fine lattice are 

ci,b = 1.068 , C2,b = 0.9234 . (24) 

For these choices we get for the coupling between the hue and the coarse lattice 


/?i,2 = 0.185807 + 0.974 (/3i - /3e) . (25) 

We recommend to use this approach in a range of couplings give by 

0.386 ^(3i^ 0.391 . (26) 

Note that 'C2n(i(0.386) = 12.5747(24). For smaller correlation length there is little need to 
use systems with varying lattice spacing. 


V. BENCHMARK: FILMS WITH STRONGLY SYMMETRY BREAKING BOUND¬ 
ARY CONDITIONS 


We intent to use the results obtained here in the study of the thermodynamic Casimir 
force between a sphere and a plane. To this end, in the neighborhood of the plane, the 
lattice spacing a = 1 is used. Also the sphere should be well inside the part of the lattice 
with lattice spacing a = 1. At increasing distances from the plane, larger and larger lattice 
spacings are used. This way, a semi-infinite system could be approximated at relatively 
low computational costs. The systematic error that is introduced by the coarsening of the 
lattice could be controlled by varying the thickness of the sectors of a given lattice spacing, 
in particular the one next to the plane. 

Here, as a benchmark and preliminary step, we study the thermodynamic Casimir force 
between two plates. This geometry is usually called film geometry. In particular we shall 
simulate systems with strongly symmetry breaking boundary conditions. This problem 
has been thoroughly studied. 


For example, in refs. [16|, [17 


Monte Carlo simulations of 

the Ising model on the simple cubic lattice were performed. In refs. 10, [l7|, llS] the im¬ 


proved Blume-Capel model on the simple cubic lattice, the model that we consider here. 
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has been studied. Hence accurate numbers to compare with are available. In addition 
to the thermodynamic Casimir force we determine the magnetization prohle. Films with 
strongly symmetry breaking boundary were also studied for example by using the extended 


de Gennes-Fisher local-functional method 19-21 


A. Thin films with strongly symmetry breaking boundary conditions 


We study a lattice of the size (Lq -|- 2) x Li x L 2 , where Li = L 2 = L throughout. 
In 1- and 2-direction periodic boundary conditions are imposed. In our convention xq G 
{0,1, ...,Lo + 1}. Strongly symmetry breaking boundary conditions are implemented by 
hxing the spins at Xq = 0 and Lq -|- 1 to either = —1 or = 1. In the case of (-F, -|-) 
boundary conditions, the spins at both boundaries Xq = 0 and Lq + 1 are hxed to = 1, 
while for (-(-, —) boundary conditions, the spins at Xq = 0 are hxed to s^, = 1 and those at 
Lq -|- 1 are hxed to s^, = — 1. In order to determine the thermodynamic Casimir force, we 
simulated two hlms of the thickness Lo,i = Lq -|- 1/2 and Lo ,2 = Lq — 1/2. In the following 
Lo ,2 is always a multiple of 4. In our simulations we used two diherent lattice spacings only. 
In particular, we kept the layers xq = 0, 1, ..., Lo,2/4 and xq = 3Lo, 2/4 -|- 1, 3Lo,2/4 -|- 2, 
... Lo,j -|- 1 as hne lattice. The remaining ones, in the center of the him, are replaced by a 
coarse lattice with twice the lattice spacing. 

The thermodynamic Casimir force per area is dehned by 

dfex 

asimir — r^ r ? 

OLq 



where fex = ffUm — LQfbuik is the reduced excess free energy per area of the him of the thick¬ 
ness Lq and ffUm is the reduced free energy per area of the him and fbuik is the reduced free 
energy density of the bulk system. On the lattice the partial derivative flTTI) is approximated 
by a hnite diherence 


dfe, 

dLr 


A/e 


(28) 


where 

A/ea;(-ho) = fex{LQ + 1/2) — fex{LQ — 1/2) . (29) 

Here, following Hucht we compute A/ea, by integrating the corresponding diherence of 
excess internal energies per area 


Mex = - / d/3AEe. 

ho 

where 

AF^ea:(-ho) = 2S.E{Lq) — Ebuik 

where E^uik is the energy density of the bulk system and 


AL(Lo) = E{Lq + 1/2) - E{Lq - 1/2) 


where we dehne the energy per area of the homogeneous system as 


E 


1 

L 1 L 2 


i^xSy) ■ 

<xy> 


(30) 

(31) 

(32) 


(33) 
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Now let us discuss how this approach can be generalized to our system with two different 
lattice spacings. The partition function takes the form 

Z= ^ exp(/3e„({s}) + ^(o. + M)e.({S},{s})--D^4--D^Si) (34) 

{S},{s} \ i X X J 

where the {S'} are the spins on the coarse lattice and {s} the ones on the fine lattice. The 
Cj correspond to sums in eq. ([2]). Taking the derivative of the reduced free energy per area 
with respect to /5 we arrive at 


df 1 d\nZ 
8(3 L 1 L 2 8/3 


L1L2 


(co) + 


which means that eq. fl55]l is generalized to 


(35) 


E 


1 

L1L2 


(co) + bi{^i) 

i 


(36) 


B. Numerical results 


Here we have repeated a part of the simulations discussed in ^ction 5 of ref. 
parameters of the simulation are very similar to those of ref. 
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The 

In particular we used 

the exchange cluster algorithm and the improved estimator that is associated with it, to 
compute AE. The implementation of the exchange cluster algorithm for the system with 
two different lattice spacings and the associated Hamiltonian is straight forward. 

Let us first discuss the results obtained for AE for (+, +) boundary conditions. We 
simulated (+, +) boundary conditions for (Lq, L) = (16.5, 64), (32.5,128), and (64.5, 256) at 
about 100 values of (3 in the neighborhood of the critical point. 

In the upper part of Fig. IHwe plot the difference 5 = AE — AE^om and 0.lAF^ex,hom for 
Lq = 16.5 and L = 64, where hom indicates the result for the homogeneous system, taken 
from ref. lOj]. We see that 6 is definitely different from zero for most of the values of (3 that 


we have simulated at. On the other hand, on average |5| is much smaller than \AEex,hom\- 
Therefore the result that we get for the thermodynamic Casimir force is at least qualitatively 
correct. Going to larger thicknesses the error due to the coarsening in the center of the film 
becomes less and less visible. While for Lq = 32.5 still |(5| 7 ^ 0 at the level of our statistical 
accuracy, for Lq = 64.5 a deviation of AE from AE^o^ can hardly be seen. See the lower 
part of Fig. |U 

Let us turn to (+, —) boundary conditions. We simulated films of the thicknesses Lq = 
16.5, 32.5 and 64.5. For each thickness we simulated at several values of L. For example 
for Lq = 16.5 we studied L = 64, 128 and 256, where L is increased with increasing f3. For 
a more detailed discussion see ref. [l^. Our numerical results for 5 = AE — AE^orn for 
Lq = 16.5 and 64.5 are shown in Fig. O For Lq = 16.5 we see a strong deviation of AE 
from the corresponding result for homogeneous system in the low temperature phase. This 
deviation is only one order of magnitude smaller than the quantity AEex,hom that we like 
to compute. For the thickness Lq = 32.5 the systematic error 5 is still clearly visible and 
its functional form is similar to that for Lq = 16.5. As for (+, +) boundary conditions, for 
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FIG. 4. We plot 6 = AE — AEhom for (+, +) boundary conditions as a function of 13. In addition we 
show as solid red line, 0.1AEex,hom- Results for the lattice size {Lq,L) = (16.5,64) and {Lq,L) = 
(64.5, 256) are given in the upper and lower part of the figure, respectively. 
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FIG. 5. Same as Fig. [3]but for (-I-, — ) instead of (-I-, -|-) boundary conditions. 
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Lq = 64.5 we find that 6 is essentially consistent with zero at the level of our statistical 
accuracy. 

As an additional check, we computed for Lq = 64.5 the integral 

A = / d(3 6 (37) 

Jpo 

which gives the error of A/. In the case of (+, +) boundary conditions, we hnd that A is 
consistent with zero within two standard deviations in the whole range of (3 that we have 
studied. The some observation holds for (+, —) boundary conditions up to = 0.391. Then 
1^1 is increasing. For {3 = 0.392 we see a deviation from zero by a little less than three 
standard deviations. 


1. Magnetization profile 


The magnetization prohle of the homogeneous system is dehned by 


m(xo) = 


L1L2 


\Xi,X2 


(38) 


Note that for (+, ±) boundary conditions 

m{xo) = ±m{Lo + 1 - xq) . 

For a discussion of the behavior at the critical point see 
results are given in section VI. B of 


18 


(39) 

and refs, therein. Numerical 
see in particular Fig. 1. For simplicity we only 
discuss the second system that we simulated, since for Lo ,2 = 16, 32 and 64 in contrast to 
Lq i = 17, 33 and 65, the symmetry fl3^ is not broken by the coarsening. 

Let us hrst discuss the results obtained for the hne part of the lattice with those for the 
homogeneous system for the same lattice size. In Fig. |6]we plot the difference 


6m{xo) = m{xo) - rrihomixo) (40) 

for Lq = 64 for both (+, +) and (+, —) boundary conditions. For individual values of (3 
we can hardly distinguish 5m{xfi) from zero. Therefore, we have averaged 5m{xQ) over all 
values 0.3852 < f3 < 0.38792 we simulated at. This way a small systematic error is revealed. 

We also determined 6m for the thicknesses Lq = 16, 32. Already here the systematical 
error is rather small. 

For the coarse part of the lattice the comparison is less direct. First one has to rescale the 
magnetization prohle by a factor of to compare it with the corresponding prohle 

of the homogeneous system. Moreover, the coarse sites are located at integer values of 
xq, while in our convention the sites of the hne lattice assume half integer values of xq. 
Therefore also an interpolation of the prohle is needed for the comparison. Here we do not 
delve into this direction but simply plot our data. In Fig. [7] we plot the magnetization 
prohle at /3 = 0.387721735 for (+, +) and (+, —) boundary conditions. We give the results 
obtained for the homogeneous system and those for the coarse lattice that we have rescaled 
as discussed above. At least by eye one can only see a deviation in the case of (+, +) 
boundary conditions at xq = 17, next to the boundary between parts of the lattice with 
diherent lattice spacings. The magnetization on the coarse lattice seems a bit too small 
compared with the one of the homogeneous lattice. 
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FIG. 6. We plot the difference 6m{xo), eq. (HH]) . between the magnetization profiles of the composite 
and the homogeneous system for both (+,+) and (+,—) boundary conditions and the thickness 
Lq = 64 at the critical temperature. We averaged 6m{xo) over all values 0.3852 < /? < 0.38792 we 
simulated at. 


VI. CONCLUSIONS AND OUTLOOK 

We aim at studying critical phenomena using lattices with adaptive lattice spacing. This 
should be useful in studies of systems with nontrivial boundary conditions. We are plan¬ 
ning to study the thermodynamic Casimir force between a planar and a spherical object 
or between two spherical objects. In critical phenomena an adaptive grid is difficult to im¬ 
plement, since the held and the couplings rescale nontrivially under a change of the length 
scale. The present work is a hrst step in this direction. We consider the universality class of 
the three-dimensional Ising model. Our starting point is the improved Blume-Capel model 
on the simple cubic lattice with lattice spacing a. The general idea is to replace the sites 
within certain sectors of the lattice by coarse sites that are separated by a lattice spacing 2a. 
This procedure might be iterated to lattice spacings 4a, 8a, ... . Furthermore these sectors 
are composed of cubes with faces that are parallel to the lattice-axis. Here we made a hrst 
step in this program. In section HTTI we worked out, how two half spaces with lattice spacing 
a and 2a, respectively, that are separated by a boundary that is perpendicular to one of 
the lattice axis, can be coupled consistently. To this end, we used hnite size scaling. The 
boundary between lattice spacing a and 2a can be viewed as a defect plane. Following refs. 
[^, [^, there is one relevant perturbation with the RG-exponent y = yt — I, where yt = 1/v, 
associated with it. Hence one would expect that it is sufficient to tune the coupling be¬ 
tween the two half spaces. However it turned out that slowly decaying corrections remain. 
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FIG. 7. We plot the magnetization profile for both (+,+) and (+,—) boundary conditions and 
the thickness Lq = 64 at the critical temperature. The data points denoted by “fine” are taken 
from simulations of the homogeneous system. The data point taken from the coarse part of the 
composite lattice are rescaled by a factor of In the case of the fine lattice, the solid lines 

should only guide the eye. 


These can be removed by tuning the couplings within the layers next to the boundary. Our 
hnal results for the coupling constants are summarized in section II V Cl In a more general 
scenario we have also edges and corners. These can be regarded as line and point defects, 
respectively. Following the argument given in the appendix of the RG-eigenvalues are 
yi = yt — 2 = —0.41275(25) and yp = yt — 3 = —1.41275(25), respectively. Hence both 
perturbations of the hxed point are irrelevant. Therefore in principle no additional tuning 
of couplings at the edges and corners is needed. However the absolute value of yi is rather 
small and hence still tuning of a coupling at the edges seems to be benehcial. It remains to 
be worked out which coupling could be tuned and which observables could be used to tune 
it. 


In order to benchmark our results we simulated dims with strongly symmetry breaking 
boundary conditions. We replaced half of the sites, located in the center of the him, by coarse 
sites. We computed the magnetization profile and the partial derivative of the internal 
energy with respect to the thickness of the him. This quantity is used to determine the 
thermodynamic Casimir force. For the largest thickness Lq = 64.5a that we have simulated, 
at the level of our statistical accuracy, the numerical results can hardly be discriminated 
from those obtained by simulating a system with a unique lattice spacing. 
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